#!/bin/bash
PANF=/home/cmbpanfs-01/yanghoch        # can be overwritten by bashrc
source ~/.bashrc
source ~/SHELL/isFileEmpty.sh
if [ $# -ne 4 ]; then
	echo "This shell output the html with 6 links to the following files"
	echo "(1) UCSC genome browser to show the coverage .bed"
	echo "(2) The candidate isoforms"
	echo "(3) The selecected isoforms"
	echo "(4) The AceView isoforms"
	echo "(5) The selected CDSs"
    echo "(6) The AceView CDSs"
	echo "Those file are expected to be generated in fixed name in the code"
	echo "The input args are:"
	echo "The 1st ptr is gene name"
	echo "The 2nd ptr is chr seq name"
	echo "The 3rd ptr is the range start"
	echo "The 4th ptr is the range end"
	echo "4 ptrs are needed"
	echo "$# is provided" 
	exit 0
fi
## Define Pathes
# (0) HTML path  
OUT_DIR=$PANF//'UCSC'//$2
OUT_LIST="GeneList$2.htm"
# (1) coverage related path
DATA_DIR=$PANF/'hg18'
WIG_PATH=$DATA_DIR/"$2.wig"             # The whole coverage info
WIG_OUT=$TMP_DIR/"$2_$3_$4.wig"         # The extracted coverage info
BED_FILE="$2_$3_$4.bed"                 # .bed file for UCSC genome browser
GENE_NAME=$(basename $1)
 
# (2)(4) selecected isoforms pathes
CANDIDATE_I=$GENE_NAME'_I.fa'
SELECTED_I=$GENE_NAME'_I_selected.fa'
SELECTED_CDS=$GENE_NAME'_CDSs_selected.fa'
# (3)(5) AceView isoforms and CDSs pathes
ACE_DIR=$PANF/AceViewMrna               # The whole AceView isoforms 
ACE_PATH=$ACE_DIR/$2mrna.fa             # The AceView isoforms in a chr
ACEVIEW_I=$GENE_NAME'_Ace_I.fa'     # Extracted AceView isoforms output path 
ACEVIEW_CDS=$GENE_NAME'_Ace_CDS.fa' # Extracted AceView isoforms output path

## Extract the files
# (1) Extract the coverage of the gene region 
#####~/SHELL/generateCoverageBED.sh $WIG_PATH $3 $4 $BED_FILE
tmp=$(python ~/SCRIPT/estimate_RPKM.py $BED_FILE 0 1000000000)
PRK=$(echo $tmp | cut -f 5 -d ' ')
known_juncs=\
$(wc -l '2_'$GENE_NAME'.junctions' | cut -f 1 -d ' ')
detected_juncs=\
$(wc -l '1_'$GENE_NAME'.bed' | cut -f 1 -d ' ')
new_juncs=\
$(wc -l 1_$GENE_NAME'.bed_'2_$GENE_NAME'.junctions_.except' | cut -f 1 -d ' ')
if [ $detected_juncs == 0 ]; then
	new_juncs=0;	
fi 
# (2)(4)
# (3)(5) Extract the AceView Isoforms
bSingleContig="F"
python ~/SCRIPT/get_seqs_from_fasta.py $ACE_PATH $1 $bSingleContig $ACEVIEW_I
if test -s $ACEVIEW_I; then
	ls $ACEVIEW_I
	echo "Call get_CDSs_from_fasta  $ACEVIEW_I $ACEVIEW_CDS"
	python ~/SCRIPT/get_CDSs_from_fasta.py $ACEVIEW_I $ACEVIEW_CDS
	#\rm $ACEVIEW_I
else
	echo "Couldn't get CDS from AceView sequence" 	
fi

## generate HTML output
if ! test -d "$OUT_DIR"; then
	mkdir -p "$OUT_DIR"
fi 
# print header 
if ! test -s "$OUT_DIR"/"$OUT_LIST"; then	
	{
		echo "<TABLE border=\"1\">"
		echo "<CAPTION><EM>Quantification of isoforms' expression levels"
		echo "<TR>"
		echo "<TH>Gene</TH>"
		echo "<TH>Expression</TH>"
		echo "<TH>RPK</TH>"
		echo "<TH>Candidate</TH>"
		echo "<TH COLSPAN=2>Selected</TH>"
		echo "<TH COLSPAN=2>AceView</TH>"
		echo "<TH>Known junctions</TH>"
		echo "<TH>Detected junctions</TH>"
		echo "<TH>New junctions</TH>"
		echo "</TR>"			
	} >> "$OUT_DIR"/"$OUT_LIST"
fi 
# print table 
URL=$(head -n 1 $BED_FILE)
URL=$(echo $URL | sed 's/^#//g')
SERVER="http://dw408k-2.cmb.usc.edu:8080/$2" 
if [ $detected_juncs > 0 ]; then
{
	echo "<TR>"
	echo "<TH>$1</TH>"
	echo "<TH><a href=\"$URL\">$2_$3_$4</a></TH>"
	echo "<TH>$PRK</TH>"
	echo "<TH><a href=\"$SERVER/$CANDIDATE_I\">Isoforms</a></TH>"
	echo "<TH><a href=\"$SERVER/$SELECTED_I\">Isoforms</a></TH>"
	echo "<TH><a href=\"$SERVER/$SELECTED_CDS\">CDS</a></TH>"
	echo "<TH><a href=\"$SERVER/$ACEVIEW_I\">Isoform</a></TH>"
	echo "<TH><a href=\"$SERVER/$ACEVIEW_CDS\">CDS</a></TH>"
	echo "<TH>$known_juncs</TH>"
	echo "<TH>$detected_juncs</TH>"
	echo "<TH>$new_juncs</TH></TR>" 
} >> "$OUT_DIR"/"$OUT_LIST"
else 
    echo "There are $detected_juncs junctions" 
fi

## move 
#cp $BED_FILE "$OUT_DIR"/"$BED_FILE"
#cp $CANDIDATE_I "$OUT_DIR"/"$CANDIDATE_I"
#cp $SELECTED_I "$OUT_DIR"/"$SELECTED_I"
#cp $SELECTED_CDS "$OUT_DIR"/"$SELECTED_CDS"
#cp $ACEVIEW_I "$OUT_DIR"/"$ACEVIEW_I"
#cp $ACEVIEW_CDS "$OUT_DIR"/"$ACEVIEW_CDS"

## generate BAT file
#BAT_OUT="$OUT_DIR"/"$1.bat"
#LOCAL_DIR="C:\\TestGenes"
#LOCAL_PATH4CDS=$LOCAL_DIR\\"$SELECTED_CDS"
#echo "start $URL" > $BAT_OUT 
#echo "perl BlatAndShow.pl human hg18 DNA query,start $LOCAL_PATH4CDS 2>$1.log" >> $BAT_OUT
exit
